MICROCOPY  RESOLUTION  TEST  CHART 

NATIONAL  BUREAU  OF  STANDAROS-1963-A 


SECURITY  CLASSIFICATION  of  This  pace  (Whan  Dote  Entarad) 

REPORT  DOCUMENTATION  PAGE 

READ  INSTRUCTIONS 

BEFORE  COMPLETING  FORM 

1.  REPORT  NUMBER  J.  GOVT  ACCESSION  NO. 

NRL  Memorandum  Report  5042  /]  fy  /^  / 

s.  recipient’s  catalog  NUMBER 

n  t 

4.  TITLE  (and  Subtltla) 

RESTRAINED  LEAST-SQUARES  PROGRAM 
(RESLSQ)  USERS  MANUAL 

S.  TYPE  OF  REPORT  A  PERIOD  COVEREO 

Interim  report  on  a  continuing 

NRL  problem 

4.  PERFORMING  ORG.  RERORT  NUMBER 

7.  author^; 

J.  L.  FI  ip  pen- Anderson,  R.  Gilardi  and  J.  H.  Konnert 

•  CONTRACT  OR  GRANT  NUMBCRfftJ 

s.  performing  organization  name  and  aooress 

Naval  Research  Laboratory 

Washington,  D.  C.  20375 

10.  PROGRAM  CLEMENT.  PROJECT.  TASK 
AREA  •  WORK  UNIT  NUMBERS 

61153N;  RR011-03-41; 
60-0303-0-3 

H.  CONTROLLING  office  name  and  aooress 

12.  report  OATS 

March  28,  1983 

is.  NUMBER  OF  PAGES 

62 

14.  MONITORING  AGENCY  name  k  ADDRESS (II  d!  Hat  ant  Item  Conlralllnt  Oltlca) 

IS.  SECURITY  CLASS,  (of  thle  report) 

UNCLASSIFIED 

1S«.  OECLASSIFICATION/OOWNGRAOING 
SCHEDULE 

I«.  DISTRIBUTION  STATEMENT  (a!  thla  Rapart} 


Approved  for  public  release;  distribution  unlimited. 


17.  DISTRIBUTION  STATEMENT  (al  Uia  atalracl  antarad  In  Black  7 0 ,  II  dlllaranl  Inca  Rapart) 


IB.  SUPPLEMENTARY  notes 


19.  KCY  WOROS  (Continue  on  rovoroo  oldm  II  n» coooory  mtd  identify  by  block  number) 

(RESLSQ) 


JO.  ABSTRACT  fCanllnua  an  tavaraa  aide  II  nacacamy  and  Idonlltr  by  black  numtae) 

,  A  restrained  least-squares  program  package  is  available  for  analysis  of  x-ray  diffraction  data  from 
«mail  and  medium  sized  molecules.  The  package  consists  of  three  basic  programs. 

The  first,  SC  ATT,  prepares  a  binary  data  file  ■which  contains  scattering  factor  information  for  each 
reflection.  The  second,  MODEL,  constructs  the  “Ideal”  molecule  towards  which  the  experimental  model 
is  refined.  The  third,  RESLSQ,  performs  the  least-squares  refinement  using  conjugant  gradient  sparse- 
matrix  methods. 

A 


DD 


FORM 
1  JAN  71 


1473 


edition  or  i  nov  os  is  obsolete 

S/B  0  102*0  14-  ABO  I 


SECURITY  CLASSIFICATION  OF  THIS  RAOE  (Whan  Data  Bnlorod) 


CONTENTS 


Introduction  . 

Program  SCATT  . 

Program  MODEL  . 

Appendices  . 

Glossary  . 

Program  RESLSQ  . 

Introduction  . 

Control  cards  (input)  . 
Subroutine  descriptions 
File  information  .  .  . 
Machine  specific  code  . 
Subroutine  arrays  .  .  . 
References  . 


1 

3 

5 

17 

26 

29 

30 
34 
40 
50 
53 
55 
60 


Accession  ?or 

NT  I  S  GRAM 
DTIC  tab 

Unannounced 

Justlfication- 


□ 


By - - - — ; - 

Distribution/ _ 

""Availability  Codes 
[Avail  and/or 


RESTRAINED  LEAST-SQUARES  PROGRAM  (RESLSQ)  USERS  MANUAL 


INTRODUCTION: 

This  refinement  package  consists  of  three  basic  programs. 

The  first  program,  "SCATT",  prepares  a  binary  reflection  file 
which  contains  scattering  factor  information  for  each  reflection. 

We  have  found  that  this  saves  time  in  the  structure  factor  calcu¬ 
lation  by  eliminating  table  look-up  and  interpolation. 

The  second,  "MODEL",  constructs  an  "ideal"  molecule  from  a 
connectivity  matrix,  covalent  radii,  and  bond  angles  supplied  by 
the  user  and  then  writes  files  of  "ideal"  interatomic  distances 
towards  which  the  "real"  model  will  be  restrained  during  refine¬ 
ment. 

The  third  program,  "RESLSQ",  performs  the  least-squares 
refinement.  It  imposes  elastic  bounds  on  all  of  the  key  geometry¬ 
defining  interatomic  distances  (Konnert,  1976;  Hendrickson  & 

Konnert,  1979)  which  were  identified  by  MODEL  (or  specified  by 
the  user)  and  also  restrains  the  thermal  factor  parameters  in 
such  a  manner  that  the  relative  motion  for  pairs  of  stereochemically 
related  atoms  does  not  exceed  physically  reasonable  values 
(Konnert  &  Hendrickson,  1980) .  Thermal  parameters  used  in 
RESLSQ  are  of  the  same  form  as  are  the  betas  in  ORFLS  (Busing, 
et  al,  1975) .  Groups  of  atoms  may  also  be  elastically  res¬ 
trained  towards  coplanarity.  How  tightly  any  particular  class 
of  restraints  is  enforced  is  at  the  user's  discretion  and  may 
be  changed  during  the  course  of  a  refinement. 

The  restrained  least-squares  procedure  uses  a  greatly 
reduced  ("sparse")  matrix;  shifts  in  parameters  are  derived  from 
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this  matrix  not  by  inversion  but  by  the  iterative  conjugate 
gradient  method  as  described  by  Konnert  (1976) .  For  a  100 
atom  problem  the  "sparse"  matrix  contains  approximately  7% 
of  the  number  of  elements  needed  for  a  "full"  r.iatrix  refine¬ 
ment.  With  increasing  numbers  of  atoms  the  size  of  the  matrix 
increases  linearly  rather  than  exponentially  as  it  does  in 
full-matrix  methods . 

The  method  has  already  proven  to  be  extremely  useful  in 
the  refinement  of  macromolecules  (Konnert,  1976;  Sielecki, 
Hendrickson,  3roughton,  Delbaere,  Brayer  &  James;  1979; 

Furey,  Wang,  Yoo  &  Sax,  1979;  Honzatko  &  Lipscomb,  1982;  Reese 
&  Lipscomb,  1982) .  The  RESLSQ  program  itself  has  been  used 
on  several  small  and  mid-sized  molecules  with  up  to  as  many 
as  213  atoms  in  the  asymmetric  unit  (Flippen-Anderson,  1982; 
Flippen-Anderson  et  al,  1982;  1983;  Karle  et  al,  1983)  and 
has  become  the  standard  refinement  program  used  in  our 
laboratory . 


PROGRAM  SCATT 


(Prepares  reflection  input  file  for  RESLSQ  -  file  20) 


Written  by 

Wayne  Hendrickson 

Laboratory  for  the  Structure  of  Matter 
Naval  Research  Laboratory 
Washington,  D.  C.  20375 


Program  SCATT  prepares  the  binary  reflection  file  read  by 
RESLSQ  (file  20) .  It  contains  scattering  factor  information  for 


each  reflection. 


CONTROL  CARDS: 

Real  Cell  Dimensions  (6F8.2) 

ACELL 
BCELL 
CCELL 
ALPHA 
BETA 
GAMMA 

Type  2  #  of  Kinds  of  Atoms  (15) 

(one  c?rd)  NS 

Scattering  factor  tables  (4x,I4,9F8.5) 

J  =  scattering  factor  order  - 

i.e.  if  J  »  1  for  carbon  scattering 
factors  then  FMULT  (FILE  10  -  atom 
coordinate  input  to  RESLSQ)  must  ■  1.0 
for  carbon. 

Cols  9-80  carry  the  9  scattering  factors 
needed  for  atom  type  J. 

Ref:  Acta  Cryst.  A24 ,  321  (1968) 
or 

Int.  Tables  -  Vol.  IV 
Table  2.2B  -  Page  99 

SCATT  also  reads  H,K,L,  FOBS  and  SIG  values  from  file  10 
(not  related  to  file  10  for  RESLSQ) 


Type  3 
(NS  cards) 


Type  1 
(one  card) 


PROGRAM  MODEL 


WRITTEN  BY: 

PAUL  B.  ANDERSON 

|  FOR: 

JUDITH  L.  FLIPPEN-ANDERSON 
LABORATORY  FOR  THE  STRUCTURE  OF  MATTER 
NAVAL  RESEARCH  LABORATORY 
WASHINGTON,  D.C.  20375 
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1.  Introduction.  The  MODEL  program  accepts  a  description  of 
the  bonds  comprising  a  molecule  and  outputs  idealized  distances 
between  various  atom  pairs . 


2.  Parsing  the  Input.  The  molecule  description  in  the  syntax 


shown  by  the  examples  of  Figure  1  is  processed  by  the  PARSE 
subroutine.  PARSE  is  repeatedly  called  by  the  BUILD  subroutine 
until  the  input  is  completely  read.  PARSE  uses  a  number  of 
supporting  routines  to  analyze  the  input  characters.  The 
characters  are  treated  as  a  continuous  stream  without  regard 
to  card  boundaries.  Routine  NXCHAR  provides  the  characters, 
one-by-one,  to  the  parsing  subroutine.  PRCHAR  is  used  in  some 
cases  to  return  characters.  PARSE  controls  the  entire  scan  with 
a  method  based  on  a  finite  state  model  of  the  input  grammar. 

The  principal  advantages  of  this  are  that  it  is  largely  table 
driven,  guaranteeing  uniformity  and  easy  maintenance,  and 
that  a  relatively  complex  structure  can  be  easily  handled.  The 
array  TRANS,  shown  in  Appendix  A,  and  the  variable  CURST  con¬ 
stitute  the  core  data  structure  of  PARSE.  CURST  gives  the 
current  "state"  of  the  parse  (see  Appendix  B) .  PARSE  analyzes 
the  input  and  determines  the  lexical  class  (see  Figure  2)  of  the 
next  characters  using  the  routines  ATOM IN  and  LISTIN.  The 
transition  table  TRANS  is  accessed  with  the  current  state  and 
the  lexical  class  of  the  input  to  determine  what  the  next  state 
will  be  and  what  action  is  to  be  taken.  Some  of  the  actions  in 
PARSE  result  in  data  being  returned  to  BUILD. 


C1:S(0.77,109.5)-H2:S(0.23,0.0) , -H3 :S , -H4 :S , -H5 : S 


represents 


H2 

I 

H4 - Cl 

I 

H5 


Cl : S-H6 : S , -H7 : S , -H8 :S 

C2:I (0.695 , 120.0) -Cl :S, -03 :S (0.66,109.5) 
C2:D (0.665, 120.0) -04:D(0.55, 120.0) 

03 :S-H5 :S 


represents 


H6  ( 


H2 


H2 


Figure  1.  Input  Syntax  Examples 


Complete  atom  eennn:b 
Incomplete  atom  sb 


1 


K- 


2 

List  (r1,r2)  3 
Comma  or  space  4 
Dash  5 
End  of  file  6 
Other  7 


Figure  2.  Lexical  Classes 


% 


i 


3.  Data  Structures  Created  by  BUILD  Phase.  Two  data  structures 
are  created  by  subroutine  BUILD.  The  first  structure  is  the 
internal  representation  of  the  bond  information  supplied  in  the 
input.  The  structure  accommodates  atoms  numbered  from  1  through 
200  and  up  to  400  bonds. 

Access  to  the  bond  information  is  through  the  atom  entries 
shown  in  Figure  3.  The  bond  entries  are  organized  as  parallel 
arrays,  storing  all  of  the  information  on  atom  "I"  in  the  Ith 
subscript  position.  Thus,  the  character  name,  bonding  angle  and 
bond  list  pointer  are  stored  in  array  elements  ATOML(I), 

AANGLE(I)  and  BONDL { I ) ,  respectively.  The  atom  number  is  the  key 
to  access  of  the  atom  list.  The  bond  list  may  be  accessed  using 
the  subscript  supplied  in  the  BONDL  entry  for  the  atom  of 
interest.  The  bond  list  is  organized  as  (five)  parallel  arrays 
in  the  same  manner  as  the  atom  list.  Figure  4  details  the 
contents  of  the  bond  list  entries.  The  bond  entries  are  chained 
together  as  a  linked  list.  The  link  between  bond  entries  is 
given  in  the  BONDNX  part  of  a  bond  entry.  The  BONDNX  value  is 
the  subscript  value  of  the  next  bond  entry  for  the  atom.  The 
list  is  terminated  by  a  zero  link.  To  speed  up  processing  of 
non-bonded  distances,  bonds  are  always  stored  twice.  If  atom 
"A"  is  bonded  to  atom  ”B",  an  entry  is  made  in  atom  A's  bond 
list  to  reflect  this.  A  corresponding  entry  is  made  in  atom  B's 
bond  list  (showing  a  bond  to  A) .  To  distinguish  the  bonds 
input  from  these  generated,  duplicate  bonds,  the  BONDNX  values 
of  the  duplicates  are  negated. 
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ATOML 


AANGLE 


BONDL 


Two  character  representation  of 
the  element  symbol,  e.g.  'CL'. 


Angle  between  bonds  for  this  atom 


Subscript  of  first  bond  entry  for 
this  atom  (the  origin  atom) . 


Figure  3.  Atom  List  Entries 


ATOMNO 


Number  of  atom  being  bonded  to 
(target  atom) . 


BNDTYl 

« 

.  ■  " 

BNDTY2 

BDIST 

rra 

J* 

v 

r  *' 

BONDNX 

m 

1  character  bond  type  for  origin 
atom. 


1  character  bond  type  for  the 
target  atom. 


Bond  distance  -  set  to  0  in  BUILD , 
computed  in  BONDD. 


Index  of  next  bond  entry  for  origin 
atom,  negated  if  this  bond  was  not 
input  (i.e.,  a  bond  created  by  BUILD). 


Figure  4.  Bond  List  Entries 


The  second  structure  built  during  subroutine  BUILD'S 
execution  stores  the  bond  distances  and  angles.  Ordinarily 
this  structure  is  retrieved  from  the  unformated  file  on  unit 
1.  User  input  may  update  or  add  to  this  data  if  desired. 

The  fina^  structure  is  printed  and  output  to  unit  2  before 
processing  completes.  To  add  to  or  update  the  structure, 
subroutine  BUILD  uses  the  auxiliary  routine  ADDRAD.  The  structure 
is  similar  to  the  atom  and  bond  structures  described  above  except 
that  the  element  list  (ELEM)  entries  are  not  ordered  in  any  way. 
Thus,  to  find  a  particular  bond  type  for  an  element,  the  element 
is  first  located  in  the  element  list.  RADL,  a  parallel  array, 
then  gives  the  index  of  the  first  bond  type  entry.  This  is 
shown  in  Figure  5.  The  bond  type  entries  are  similar  in 
organization  to  the  bond  list  entries  of  the  bond  data  structure. 
Bond  type  entries  have  the  structures  shown  in  Figure  6.  (Note: 
Structure  here  is  used  in  a  programming  sense  and  does  not  refer 
to  the  "model  structure"  being  described  by  the  restraint  para¬ 
meters.  ) 

4.  Bonded  Distances.  The  bond  distances  are  computed  by 
routine  BONDD  once  the  input  has  been  completely  scanned. 

BONDD  also  establishes  the  bond  angles  used  in  subroutine 
NBONDD  to  compute  the  non-bonded  distances.  Bond  distances  are 
determined  by  the  bonding  radii  of  the  two  atoms  involved. 

These  bonding  radii  vary  for  each  element  depending  on  the  bond 
type  (e.g.  single  or  double).  The  bond  distance  is  computed 
simply  as  the  sum  of  the  two  radii.  The  bond  distances  are 
output  to  unit  26  in  the  format  shown  in  Figure  7.  Each  bond 
distance  is  stored  in  the  bond  entry  (BDIST)  for  use  by  NBONDD. 


ELEM 


RADL 


2  character  symbol  for  element, 
e.g.  'CL'. 


Subscript  of  first  bond  type  entry 
for  this  element. 


Figure  5.  Element  List 


BTYPE 


1  character  bond  type. 


BONDR 


ANGLE 


BRLNXT 


Radius  for  a  bonded  distance 
of  this  type. 


Angle  for  bonds  (used  to  set 
array  AkNGLE) . 


Subscript  of  next  bond  type 
entry.  0  if  no  more  entries 


Figure  6.  Bond  type  Entries. 


Each  record  is  produced  by  unformatted  I/O  in  the  following  format: 


P  4J 
0  0)  — . 
4J  O'  ^ 
(d  P  0) 


P  <d  O' 
0)  EH  <D 

a  -p 


© 

in 


O 

•P 

■P 

a> 

(0 

a> 

P 

(0 

in 

C 

o 

-H 

■P 

id 

rH 

CQ 

c 

id 

p 

+> 

d) 

o 

-P 

•P 

•P 

id 


0) 

<u 

• 

O' 

a) 

rH 

a) 

p 

■p 

A 

O 

c 

H 

4J 

■p 

W* 

rH 

-p 

rH 

a) 

P 

a) 

< 

w 

■fl 

s 

+ 

* 

* 

8 

■p 

< 

Z 

H 

a 

H 

PS 

o 


0) 
r— I 

3 

O 

a> 

iH 

O 

s 

>1 

-P 

-H 

■P 

c 

a) 

•o 


a) 

a 

4J 

P 

O 

IP 

+» 

id 

e 

p 

o 

ip 

-p 

3 

a 

•p 

3 

O 

ID 

IN 

•P 

-H 

c 

p 


<u 

p 

3 

O' 

•H 

IP 
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5.  Non-Bonded  Distances.  Non-bonded  distances  are  computed 
using  the  law  of  cosines  for  any  two  distinct  atoms  bonded  to 

a  common  atom.  The  separation  of  these  atom  pairs  is  determined 
by  the  bond  distances  and  the  bond  angles  created  by  BONDD. 
Non-bonded  distances  are  output  to  unit  26  in  the  format  given 
in  Figure  7 . 

6.  The  user  may  supply  specific  distance  restraint  information. 
Distance  information  would  be  input  through  unit  16  and  added 
by  SUBROUTINE  PREFIX  to  unit  26  output. 

7.  Block  Data  Output.  The  bonding  radii  and  angles  used  for 
BONDD  and  NBONDD  are  kept  in  common/RADII/.  This  common  is 
created  by  the  user  input  and/or  read  from  unit  1  when  the 
program  first  begins  execution.  Once  the  program  has  completed 
the  output  for  unit  26  as  described  above,  the  RADII  common 

is  output  using  unformated  I/O  to  unit  2  so  that  it  may  be 
retained  for  subsequent  use  (as  unit  1) .  The  bond  descriptions 
are  thus  easily  expandable  and  can  retain  their  values  across 
program  executions.  It  is  intended  that  the  tables  in  RADII 
reach  a  point  where  additions  are  few  and  those  that  do  occur 
are  for  special  situations. 

APPENDIX  D  contains  file  descriptor  information,  details 
on  control  cards  necessary  to  execute  the  program  and  a  glossary 
of  terms  which  are  used  in  the  COMMON  blocks. 
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Appendix  A 
TRANS (*»*,!) 


TRANSITIONS 


1  2 

1  2  0 

2  0  2 

3  4  0 

4  0  0 

5  2  2 

6  6  6 

7 

8 
9 

10 


3  4  5  6 
0  10  0 
2  0  3  0 
0  0  0  0 

4  5  0  6 
0  5  3  6 
6  6  6  6 


7  8  9  10 

0  0  0  0 

0 
0 
0 
0 
6 

0 

0 

0 

0 


1 

2 

3 

4 

5 

6 

7 

8 
9 

10 


TRANS (*, 

*,2) 

- 

ACTIONS 

1 

2 

3 

4 

5 

6 

7  8 

1 

0 

0 

7 

0 

0 

0  0 

0 

4 

2 

0 

3 

0 

0 

5 

0 

0 

0 

0 

0 

0 

0 

0 

6 

7 

0 

11 

0 

9 

10 

0 

7 

8 

11 

0 

12 

12 

12 

12 

12 

12 

12 

0 

0 


0 

0 


NOTE:  Vector  elements  which  have  been  left  blank  are  zero. 
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APPENDIX  B 


Appendix  C 


MACHINE  DEPENDENT  ROUTINES 

The  following  routines  may  require  modification  for 
other  machines  than  the  TI  ASC.  In  addition,  any  routines 
performing  I/O  may  be  susceptible  to  change. 

NXCHAR 

NUMER 

ALPHA 

PACK2 

VALUE 

BUILD 


APPENDIX  d 


FILE  INFORMATION 


FILE  (01) 


FILE  (02) 


File  (05) 
File  (06) 
FILE  (16) 
FILE  (26) 


Contains  a  dictionary  of  information  on  distances 
and  angles  to  generate  the  model.  (E.g.  an  entry 
of  C1:S (0 . 77 , 109 . 5)  would  create  an  entry  for  an 
atom  C  with  bond  type  S  having  a  covalent  radius  of 
0.77A  and  a  covalent  angle  of  109.5).  This 
dictionary  can  be  added  to  or  changed  any  time  the 
program  is  run. 

Contains  an  updated  verstion  of  the  dictionary  in 
file  01  and  is  what  is  output  to  save  as  the  new 
file  01  for  a  future  run. 

If  file  manipulation  is  a  problem,  these  files  can 
be  eliminated  and  the  dictionary  information  for 
each  model  can  be  supplied  as  part  of  the  input 
date. 

Control  card  input. 

Printout. 

User  supplied  distance  constraint  information. 

Output  file  for  distance  constraint  parameters. 

(See  Figure  7.) 


INPUT  INFORMATION 


i 


FILE  (01)  File  01  is  a  binary  file  which  has  been  written 

(to  file  02)  and  saved  from  an  earlier  MODEL  RUN. 

FILE  (16)  FORMAT ( 215 ,F7 . 3 , 3X , II , 4X , 3F10 . 4 , 3X , 12) 

ANOl  -  First  atom  number 

AN02  -  Second  atom  number 

DIST  -  Distance  between  ANOl  and  AN02 

TYPE  -  Integer  code  for  distance  type 

1  »  Bonded  (non-hydrogen) 

2  =  Non-bonded  (non-hydrogen) 

3  =  Bonded  (involving  a  hydrogen) 

4  =  Non-bonded  (involving  a  hydrogen) 

5  =  Intermolecular  distance 

XSYMM1 

XSYMM2  -  lattice  translations  for  target  atom 
XSYMM3 

LSYMl  -  Symmetry  operator  number  for  target  atom 
NOTE:  For  user  added  distances  it  will  be  necessary 

to  set  the  location  of  the  "target"  atom  by  defining 
a  set  of  lattice  translations  (XSYMM1,  XSYMM2 ,  XSYMM3) . 
For  a  lattice  translation  of  one  unit  cell  along  a, 
XSYMMl  will  be  6.0,  etc.  A  symmetry  operator  number 
(LSYMl)  must  also  be  defined.  This  is  the  numerical 
position  of  the  symmetry  operation  in  the  input  stream 
for  RESLSQ .  These  operations  will  be  carried  on  the 
coordinates  of  the  target  atom  as  it  appears  in  the 
input  stream. 

FILE  (05)  CONNECTIVITY  DESCRIPTOR  CARDS  (FREE  FORMAT' 

**  See  additional  description  on  following  pages. 


21 


OUTPUT  INFORMATION 


FILE  (02)  Dictionary  of  bonding  radii  and  angles  -  updated 
and  output  after  processing  completed. 

FILE  (06)  Printout. 

FILE  (26)  Distance  constraint  information. 


FILE (05)  CONNECTIVITY  DESCRIPTOR  CARDS  (FREE  FORMAT) 


e.g.  C1:S(0.77,109.5)-C2:S,-H3:S(0.23,0.0) ,-H4:S 

C  =  atom  name  -  2  characters  are  allowed  and  any  name 
the  user  wishes  is  allowable  as  long  as  one  is 
consistent.  Hydrogen  atoms,  however,  must  be  referred 
to  as  H  because  the  program  keys  on  that  character  to 
flag  the  different  bond  types. 

1  =  position  atom  will  occupy  in  CONLSO  input  deck.  Any 
numbering  scheme  is  allowable  as  long  as  it  is  con¬ 
sistent  with  atom  input  deck  to  CONLSQ.  A  number  may 
be  used  only  once. 

S  =  1  character  code  for  type  of  bond  associated  with  a 
particular  atom  name  (26  different  bonding  environ¬ 
ments  are  possible  for  each  atom  name) . 

(0.77,109.5)  -  the  bonding  radius  and  bonding  angle  for  the 
particular  atom  and  bond  type  immediately  preceding  this 
information.  Need  be  put  in  only  once  for  each  different 
atom  and  type  of  bond.  A  dictionary  of  these  values  can 
be  built  and  used  in  later  runs. 


DELINEATORS 


COLON  (:) 


PARENTHESES  ( ) 


HYPHEN  (-) 


COMMA  ( , ) 


BLANK 


Separates  an  atom  name  and  number  from 
a  bond  type  indicator. 

Enclose  dictionary  information.  This 
information  should  be  supplied  from  the 
05  input  at  the  first  occurrence  of  the 
atom-bond  type  it  describes  (or  from  unit  01) 
Separates  first  and  second  atoms  describing 
a  bond. 

Indicates  information  about  another  bond 
follows 

A  blank  card  ends  the  input  deck. 


NOTE:  It  is  only  necessary  to  describe  (in  the  manner 

illustrated  above)  each  bonded  distance  in  the 
molecule.  MODEL  will  determine  non-bonded  re¬ 
lationships  and  set  thermal  parameter  constraints 
(according  to  Figure  9) . 


SAMPLE  INPUT  DECK 


C1:S(0.77,109.5)-C2:S,-09:S(0.66,109.5) , -012 : S , -H16 : S ( 0 . 23 , 0. 0) 

C2:S-C3:S,-N15:S(0.70,109.5) ,-H17:S 

C3:S-C4:S,-OlO:S,-H18:S 

C4 : S-C5 : S , -Oil : S , -H19 : S 

C5 : S-C6 : S , -012 : S , -H2  0 : S 

C6 : S-013 : S , -H21 :S , -H22 : S 

C7:I(0. 695,120.0) -N15 : S , -C8 : S 

C7 : D (0. 665 , 120. 0) -014 :D (0.55 ,120.0) 

C8 : S-H2 3 : S , -H2 4 : S , -H25 : S 


09 :S-H26 :S 
OlO:S-H27:S 
Oll:S-H28:S 
013 :S-H29 : S 
N15 : S-H30 : S 


(R.  D.  Gilardi  &  J.  L.  FlipDen  (1974), 
Acta  Crvst . ,  B30,  2931) 
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GLOSSARY 


COMMON  ALIST 

ATOML  -  (Atom  list)  2-character  element  names  for  atoms 
in  use.  Set  in  BUILD. 

BONDL  -  (Bond  list)  subscript  of  first  bond  entry. 

Set  in  BUILD. 

AANGLE  -  (Atom  Angle)  angle  of  bonds  for  a  given  atom. 
Set  in  BONDD. 


COMMON  BONDS 
ATOMNO 

BNDTY1 

BNDTY2 

BDIST 

BONDNX 

SYMMl 
SYMM2 
SYMM3 

LSYMM 


(Atom  number)  atom  numbers  of  target  atom  in 
bond  entry.  Set  in  BUILD. 

(Bond  type  1)  Bond  type  for  originating  atom. 

Set  in  BUILD. 

(Bond  type  2)  Bond  type  for  target  atom. 

Set  in  BUILD. 

(Bond  distance)  bond  distance.  Set  in  BONDD. 

(Bond  next  )  subscript  of  next  bond  entry. 

Set  in  BUILD. 

Sets  position  of  identity  molecule  relative  to 
translations  along  cell  edges.  All  three  are  set  to 
5.0  which  is  equivalent  to  a  lattice  translation  of 
0.0  along  all  cell  edges  for  the  origin  molecule. 

=  1  (symmetry  operator  #  for  origin  molecule) 


COMMON  INPUT  -  Used  primarily  in  routines  PARSE  and  BUILD. 

ANO  -  (atom  number  )  Atom  number  of  originating  atom 

ANM  -  (atom  name)  2-character  element  symbol. 

-  (Bond  type  1)  1-character  bond  type  for  ANO. 


BT1 


BNO 

BNM 

BT2 

NBNM 

NBT 

NBDIST  - 
NBANG 
CHAR 
LPTR 

LINE 

COMMON  RADII 
ELEM 
RADL 

BTYPE 
BONDR 
ANGLE 
BRLNXT  - 

COMMON  SPS 
BLSP 


(Bonding  atom  number)  atom  number  of  target  atom. 
(Bonding  atom  name)  2-character  element  symbol 
of  target  atom. 

(Bond  type  2)  1-character  bond  type  for  BNO. 

(new  bond  name)  2-character  symbol  for  new  bond 
types . 

(new  bond  type)  1-character  bond  type  for  new  bond 
types. 

(new  bond  distance)  bond  radius  for  new  bond  types, 
(new  bond  angle)  angle  for  new  bond  types, 
(character)  current  input  character. 

(line  pointer)  current  position  (column)  in  input 
line. 

(current  input  line)  current  input  line. 


(Element  name)  2-character  element  symbol. 
(Radius  list)  subscript  of  first  bond  type  entry 
for  each  entry  in  ELEM. 

(bond  type)  1-character  bond  type 
(bonding  radius)  bonding  radius 
(bonding  angle)  bonding  angle 

(bond  radius  list  next)  link  to  next  bond  radius 
list  entry. 


(bond  list  stack  pointer)  Allocation  pointer  for 
bond  list  entries. 
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ASP 


-  (Atom  list  stack  pointer)  Allocation  pointer 
£or  element  list. 

RADSP  -  (Radius  list  stack  pointer)  Allocation  pointer 
for  bonding  type  list. 


COMMON  STATE  -  Used  in  PARSE 
CURST  -  (current  state) 
TRANS  -  (transition  matrix) 


COMMON  COUNTS 
BDCNT 
NBDCNT 
XBONDS 
XNBND 
XINTER 


-  Used  in  BONDD ,  NBONDD ,  and  PREFIX 

-  Number  of  bonded  distances  (No  H) 

-  Number  of  non-bonded  distances  (No  H) 

-  Number  of  bonded  distances  (with  H) 

-  Number  of  non-bonded  distances  (with  H) 

-  Number  of  intermolecular  distances 
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PROGRAM  RESLSQ 


A  restrained  least-squates  refinement  program  using 
conjugant  gradient  sparse-matrix  methods. 


Written  by 

John  H.  Konnert 
Richard  Gilardi 
Judith  L.  Flippen-Anderson 

Laboratory  for  the  Structure  of  Matter 
Naval  Research  Laboratory 
Washington,  D.  C.  20375 


RESLSQ 


INTRODUCTION 

A  least-squares  refinement  which  imposes  restraints  on 
atomic  positions  and  on  thermal  parameters  has  proven  to  be 
extremely  useful  for  macromolecules  (Konnert,  Acta  Cryst  (1976) 
A32,  614;  Hendrickson  and  Konnert,  "Aspects  of  Biomolecules" 
(1978)  Ed.  by  R.  Srinvansan,  Vol.  1,  Pergamon  Press;  Konnert 
and  Hendrickson,  Acta  Cryst  (1980)  A36) .  The  restraints 
effectively  reduce  the  number  of  parameters  to  be  determined 
which  means,  for  macromolecules,  that  meaningful  refinements 
can  be  carried  out  with  very  limited  data  sets.  Intermediate 
sized  systems  (up  to  several  hundred  atoms) ,  with  high  resolution 
data  sets,  generally  refine  smoothly  to  a  definite  endpoint  with 
full-matrix  least-squares  techniques.  However,  full  matrix 
refinement  with  anisotropic  thermal  motion  is  a  large  and 
expensive  procedure  for  mid-sized  molecules.  The  development 
of  RESLSQ,  an  anisotropic  sparse-matrix  refinement  program  with 
distance,  angle,  and  thermal  variance  restraints,  provides  a 
rapid,  economical  alternative. 

Rapid  convergence  is  due,  for  the  most  part,  to  the  fact 
that  only  elements  of  the  derivative  matrix  related  to  the 
physical  restraints  need  be  calculated  and  stored  (on  the  order 
of  1%  of  the  number  of  terms  needed  for  a  full-matrix  refine¬ 
ment  for  a  100  atom  problem) .  Storage  requirements  increase 
linearly,  not  quadratically ,  as  the  number  of  atoms  being  refined 
.  In  RESLSQ,  all  parameters  are  refined  simultaneously, 


increases 


but  only  correlations  between  bonded  atoms  and  between  atoms 
bonded  to  a  third  atom  are  normally  included.  Since  the  off- 
diagonal  elements  which  represent  these  "bonding  correlations" 
are  dispersed  (sparsely)  throughout  the  matrix,  the  problem 
is  not  separable  into  a  number  of  smaller  sets  of  equations,  as 
it  is  in  other  sparse-matrix  techniques  such  as  diagonal  or 
block-diagonal  least-squares  refinement.  Matrix  inversion  is 
not  a  space-efficient  technique  for  sparse  matrices;  the  inverse 
of  a  sparse  matrix  is  not  sparse.  However,  the  conjugate 
gradient  technique  is  an  iterative  technique  which  retains 
the  sparse  derivative  matrix  throughout  the  solution  of  the  normal 
equations.  In  the  programming,  only  the  non-zero  elements  of 
the  (hypothetical)  NxN  coefficient  matrix  are  packed  into  a  one¬ 
dimensional  program  array.  A  simple  "bookkeeping"  scheme 
enables  the  location  and  retrieval  of  any  element  as  needed. 

Hydrogen  atoms  may  be  included  with  anisotropic  thermal 
factors  from  the  onset  of  the  refinement,  since  the  thermal 
factor  restraints  ensure  that  the  thermal  parameters  of  the  model 
will  stay  within  reasonable  stereochemical  limits. 

In  our  laboratory  with  mid-sized  (approx.  100  atoms)  aniso¬ 
tropic  molecules,  we  have  achieved  an  increase  in  refinement 
cycle  speed  of  up  to  30  times  that  of  full  matrix  techiques. 
(Although  RESLSQ  may  take  more  cycles  to  converge  than  full-matrix 
methods)  We  have  also  been  able  to  do  simultaneous  anisotropic 
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refinement  on  all  atoms  in  cases  where  the  full-matrix  would 
be  too  large  to  store  on  even  a  very  large  computer. 

This  refinement  package  consists  of  two  basic  programs. 

The  first,  "MODEL",  constructs  an  "ideal"  molecule  from  a 
connectivity  matrix,  covalent  radii,  and  bond  angles  supplied 
by  the  user  and  then  writes  files  of  "ideal"  interatomic 
distances  towards  which  the  "real"  model  will  be  constrained 
during  refinement. 

The  second  program,  "RESLSQ" ,  performs  the  least-squares 
refinement.  It  imposes  elastic  bounds  on  all  of  the  key  geometry 
defining  interatomic  distances  which  were  identified  by  MODEL 
(or  specified  by  the  user)  and  also  restrains  the  thermal  factor 
parameters  in  such  a  menner  that  the  root  mean  squares 
amplitudes  of  relative  motion  for  these  atom  pairs  do  not 
exceed  physically  reasonable  values;  i.e.  0.05A  for  covalently 
bonded  atoms  and  0.1A  for  next  nearest  neighbor  atoms.  How 
tightly  these  restraints  are  enforced  is  at  the  user's  discretion 
and  may  be  changed  during  the  course  of  a  refinement. 

PHYSICAL  RESTRAINTS  IMPOSED  ON  THE  ATOMIC  MODEL: 

1.  Covalent  bonds  are  distributed  about  "ideal"  values. 

2.  Bond  angles  remain  near  "ideal"  values. 

3.  Approximately  rigid  (planar,  etc.)  groups  remain  so. 

Both  positional  and  thermal  parameters  must  be  consistent 

with  the  above  conditions. 


THERMAL  RESTRAINTS  IMPOSED  ON  THE  ATOMIC  MODEL: 


1.  The  variance  of  interatomic  distance  distributions  that 
are  related  to  bond  lengths  and  angles  are  restrained  to  suitable 
low  values.  These  variances  are  expressed  in  terms  of  individual 
atom  thermal  parameters. 

NOTE:  Throughout  the  program  we  are  dimensioned  for  500 

atoms  -  except  for  AM,  AMM,  A  and  S  (the  latter  are 
in  SCALE  and  SCALEB)  which  are  currently  dimensioned 
(40000)  large  enough  for  -300  atoms  and  -900  distances. 

For  more  atoms  or  distances  these  arrays  must  be  increased. 
They  require  (21  x  NA)  +  ( 36  x  NDIS)  +  1 
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RESLSQ  CONTROL  CARDS 


Type  1 : 


Type  2 : 


1-9 
10-18 
19-27 
28-36 
37-45 
46-54 
Type  3 : 
1-5 
6-10 


Type  4 : 


lations 


(one  card)  FORMAT  18A4 

TITLE  CARD  -  Any  alphanumeric  information 
(one  card)  FORMAT  (8F9.4) 

Real  cell  dimensions 
A 
B 
C 

ALPHA 

BETA 

GAMMA 

(one  card)  FORMAT  (215, 5A4) 

NSYM  =  number  of  symmetry  operations 
NCENTR  =  0  for  centro  space  groups 

1  for  noncentric  space  groups 
SPCGRP  =  20  character  field  for  comment  on  space 
group  (used  only  for  printout) 

(three  cards  per  symmetry  operation)  FORMAT  (4F5.0) 
Card  one  x 

Card  two  y 

Card  three  z 

Each  card  contains  an  x,y,z,  and  translation  (trans- 
are  in  12th' s) . 


e.g.  X  1/2  +  y  z 
Card  one  -1000 
Card  two  0106 
Card  three  00-10 


NOTE:  If  NCENTR=0 ,  omit 

symmetry  operations 
related  by  inversion 


Type  5 

(one  card) 

FORMAT  (815) 

1-5 

NA 

= 

number  of  atoms 

6-10 

NTYP 

= 

number  of  types  of  atoms 

11-15 

NDIS 

s 

total  number  of  distances  to  be  restrained 

16-20 

NPLN 

= 

number  of  planes 

21-25 

NV 

= 

number  of  variables  (NA*3  +  number  of  scale 
factors) 

26-30 

NOBS 

= 

maximum  number  of  reflections  to  be  used  in 
refinement 

31-35 

NSCAT 

s 

scattering  factor  table  number  for  H  atoms 

36-40 

IANOM 

= 

flag  to  read  anomalous  dispersion  corrections 

0  -  skip 
1  -  read 


Type  5A:  (two  cards)  FORMAT  (8F5.3) 

Read  these  cards  only  if  IANOM  =  0 
First  card  -  f'  for  each  atomic  scattering  type 
Second  card  -  f"  for  each  atomic  scattering  type 
Type  6:  (one  card)  FORMAT  (1015) 

1-5  NC  =  number  of  cycles  to  be  done  with  conditions 

currently  being  described 

6-10  NOUT  =  an  integer  value  such  that  reflections  having 

(Yo-Yc)  greater  than  NOUT  times  SIGMAF  will 
be  considered  outliers  and  will  not  be  included 
in  the  refinement 

11-15  NSK  =  0,  add  products  of  structure  factor  derivatives 

to  the  distance  blocks 


1 ,  do  not  add . . . 
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16-20 


IPRF 


21-25  IORFE 


26-30  IWGT 


31-35  IDIFF 


36-40  NRFAC 


41-45  I DUMP 


46-50  IHRID 


0,  do  not  print  reflections 

1,  print  reflections 

2,  print  only  outliers 

punch  options  for  conversion  of  Lhermal  factors 
to  form  suitable  for  ORFFE  or  for  journal  publication 
0,  no  conversion 

1,  punch  a  deck  of  S.D.'s  suitable  for  input  to  ORFFE 
-1,  product  a  listing  of  coordinates  and  big  Bij's 
as  well  as  ORFFE  deck  (coordinates  as  used  through¬ 
out  are  compatible  with  ORTEP) 

0  for  equal  weights  (i.e.  all  sigmas  on  file  20 
set  equal  to  SIGMAF) 

1  for  weighted  refinement  (i.e.  all  sigmas  on 
file  20  multiplied  by  SIGMAF) 

0,  skip  writing  difference  file 

1,  write  file  for  difference  map  on  file  30 
excluding  outliers 

2,  write  file  for  difference  map  on  file  30 
including  outliers 

0,  refine  as  usual 

1,  calculate  R-factor  only  -  no  refinement 
0,  shortened  printout 

1,  gives  dump- type  printout  of  first  few  matrix 
elements  (for  both  distt’joes  and  thermal  para¬ 
meters)  after  restraints  have  been  added 
0,  refine  as  usual 

1,  set  H  atoms  B's  to  be  the  same  as  those  of 
the  atoms  to  which  they  are  bonded 
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57-64  WPSKAL  =  overall  multiplier  for  planar  restraint 

contributions 


65-72 

SIGP 

"ideal" 

standard  deviations  from  planes 

73-80 

WBSKAL 

overall 

multiplier  for  thermal  parameter 

restraints 

Type  9: 

(one  card) 

FORMAT 

(4F8.3) 

SIGB(N)  = 

"ideal" 

motion 

root  mean  square  amplitude  of  relative 

(by  distance  classes  N-as  defined  above) 

Type  10: 

(NPLN  cards) 

FORMAT  (1615) 

1-5 

NOPLAN 

number 

of  plane  to  be  described 

6-10 

NAP 

number 

of  aroms  in  this  plane 

NAT (NAP)  = 

ate’-'  numbers  of  atoms  in  this  plane  (acorn 

number  refers  to  the  position  of  the  atom  in 
the  input  coordinate  file) 

NOTE:  Current  dimensions  allow  a  maximum  of  30  planes  with  up  to 

20  atoms  per  plane 

Type  11:  (three  cards)  FORMAT  (10F5.3) 

damping  factors  and  weighting  factor  for  some  distance 
constraints 

1-5  Card  one:  DAMPSC  =  damping  factor  for  scale  factor 

6-10  WTH  =  weight  to  be  given  to  non-hydrogen 

distance  restraints  in  the  vector 
set  up  (should  be  between  0.0  and  1.0). 

A  value  of  0.0  essentially  removes  the 
distance  restraints  -  lets  bonds  lengths 
go  where  structure  factors  take  them. 
However,  restraint  information  still 


goes  into  matrix  setup  which  effectively 
keeps  the  bond  distances  elastically 
bound  to  their  current  values  rather 
than  to  a  target  distance.  (Tends  to 
refine  by  changing  torsions,  keeping  dis¬ 
tances  and  angles  somewhat  fixed) 

NOTE:  If  DAMPS C  is  negative  then  it  will  be  applied  as  an  overall 

damping  factor  and  the  damping  parameters  on  the  following 
two  cards  will  be  ignored. 

Card  two:  DAMPC(N)  =  damping  factors  for  coordinate  shifts 

N  is  the  scattering  factor  reference 
number  (i.e.  a  different  damping  is 
allowed  for  each  atom  type) 

Card  two:  DAMPT(N)  =  damping  factor  for  thermal  parameter 

shifts  (N  has  same  meaning  as  above) 
These  cards  must  be  included  even  if  blank 

NOTE:  It  will  probably  be  necessary  to  apply  some  damping  through¬ 

out  the  course  of  the  refinement  (especially  for  thermal 
parameters)  due  to  the  high  correlation  between  scale  factor 
shifts  and  thermal  factor  shifts. 

Type  12:  Continuation  or  termination  card 
Read  another  type  6  card  here. 

If  NC=0  (or  blank)  -  program  stops 

If  NC^fO  then  a  new  set  of  refinement  conditions  will  be 
accepted  and  this  card  must  be  followed  by  type  7,  type  8, 
type  9,  and  type  11  cards  and  a  new  type  12  card  (it  is 
not  necessary  to  repeat  the  plane  information) 
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BCNVRT 


LIST  OF  SUBROUTINES 


BREF 

CALC 


CART 


CGSOLV 

CGSOL1 

COMBB 

RESLSQ  (MAIN) 

DIAG 

DISRF 

ELAPSD 

FLAG 

FREEZE 

HRIDE 

INIT 

MATMLL 

MATMUL 


PLANE 


PLNRF 

RECYCL  (ENTRY  IN  INIT) 
REDIST  (ENTRY  IN  DISRF) 
RFAC 

RSCAL  (ENTRY  IN  SCALE) 
SCALE 


SCALE  B 


SYMMIN 

TEST 
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COMMON  BLOCKS 


ATOMS 

FILES 

RECIP 

TYMNG 

CONST 

HRIDE 

VAR 

CONTROL 

LSQ1 

REFLCT 

VECTOR 

DERIVB 

LSQMTX 

RFAC 

WEIGHT 

DERI VS 

OTHER 

STDDEV 

DIFCUT 

PLANES 

SYM 

DIST 

REAL 

TABLES 

1.  BCNVRT  is  used  to  calculate  thermal  factors  which  are  refer¬ 
red  to  conventional  crystal  axes,  as  opposed  to  RESLSQ's  B's. 

In  RESLSQ  the  refinement  is  carried  out  with  thermal  factors 
(B's)  which  fit  the  form. 

exp  [-B11h2+B22k2+B33l2+2B12hk+2B13hl+2B23kl)] 


which  are  the  same  as  ORFLS  B's.  If  I0RFE=-1,  thermal  factors 
are  converted  to  the  form 


exp  [-l/4(B11h  a  +B22k 


2  *2 


+B33l2c*2+2B12hka*b*+2B13hla*c*+2B23klb*c*) 1 


which  we  sometimes  call  "big  Bij's".  These  are  the  coefficients 
currently  requested  by  Acta  Cryst.  They  are  on  the  same  scale  as 
the  usual  isotropic  B's.  BCNVRT  also  produces  a  listing  of  coordinates 
and  BEQ  values  suitable  for  submission  with  a  publication.  If  I0RFE=1, 
a  file  of  standard  deviations  suitable  for  use  with  ORFFE  is  produced. 
ORFFE  is  the  Oak  Ridge  program  for  calculating  distances,  angles, 
torsions,  etc. 
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2.  BREF  calculates  the  variance  in  each  of  the  restrained 
interatomic  distances  due  to  the  difference  in  thermal  motion 
indicated  by  the  thermal  ellipsoids  of  the  two  atoms.  (Details 
of  this  calculation  are  published:  Konnert,  J.  and  Hendrickson, 

W.  A.,  Acta  Cryst.  A36  (1980)  344-350. 

Each  of  these  variances  is  weighted  in  accordance  with 
the  expected  variance  of  its  "distance  class",  and  added  to  the 
overall  residual  being  minimized  in  the  refinement. 

The  variance  restraint  is  a  physically  reasonable  restriction, 
based  upon  experience  with  measurements  of  bond  variances  (e.g., 
electron  diffraction) .  The  variance  restraint  could  be  used  with 
any  thermal  expression.  BREF  prints  out  of  a  summary  of  actual 
and  target  variances,  with  large  discrepancies  flagged. 

3.  CALC  is  basically  a  structure  factor  calculation  routine. 

It  is  entered  once  for  each  reflection.  It  is  called  by  SF.  CALC 
uses  table  look-ups  to  get  exponentials  and  sines  and  cosines. 

This  will  probably  be  quicker  on  most  computers  than  using  the 
actual  functions  each  time. 

4.  CART (XX,X,NTC)  takes  NTC  3-component  vectors,  X,  expressed 
in  unit  cell  fractional  components,  and  transforms  them  to  NTC 
vectors,  XX,  which  have  Cartesian  Angstrom  coordinate  components. 

5.  CGSOLV  solves  the  set  of  linear  equations,  known  as  the 
normal  equations,  which  are  obtained  by  expanding  the  overall 
residual  in  a  Taylor  series  (linearizing) ,  and  minimizing  with 
respect  to  shifts  in  the  parameters.  In  matrix  shorthand,  this 
set  is  written 
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In  the  program,  these  letters  are  used  to  name  the  corres¬ 
ponding  arrays:  X(3500)  is  the  shift  vector,  R(3500)  is  the 
righthand-side  vector,  and  A (30, 000)  is  the  coefficient  matrix  - 
a  sparse  2-dimensional  matrix  whose  non-zero  elements  are  stored 
in  a  one-dimensional  computer  array.  The  normal  equations  are 
solved  by  the  conjugate  gradient  procedure  (Konnert,  J. ,  Acta 
Cryst.  (1976)  A32,  614-617;  Hestenes,  M.  R.  and  Stiefel,  E., 
Journ.  Batl.  Bur.  Stand.  (1952)  £9,  409-436.  The  notation  used 
in  CGSOLV  exactly  parallels  that  used  in  the  cited  references. 

The  conjugate  gradient  iterative  procedure  has  been  found 
to  converge  more  rapidly  if  the  normal  equations  are  scaled  in 
such  a  way  that  the  diagonal  elements  are  equal.  This  is  auto¬ 
matically  done  before  starting  the  c.g.  iterations  by  calling 
subroutine  scale.  The  unsealed  shift  is  retrieved  after  the  c.g. 
procedure  by  calling  RSCAL  (an  alternate  entry  in  SCALE) . 

CGSOLV  is  "special"  rather  than  general  mainly  because  of 
the  bookkeeping  scheme  necessary  when  storing  just  a  few  2- 
dimensional  components  in  a  1-dimensional  array.  It  will  have 
to  be  modified  whenever  this  scheme  is  modified  (e.g.,  if  more 
or  less  thermal  parameters  are  used,  or  more  scale  factors,  etc.) 
6.  CGSOL1  is  almost  exactly  the  same  as  CGSOLV,  and  one  or 
the  other  could  probably  be  dispensed  with  by  recoding.  The 
thermal  parameter  refinement  is  completely  separable  from  the 
position  refinement  and  certain  subroutines  (CGSOL1,  MATMLL) 
deal  only  with  the  thermal  refinement. 


7.  COMBB  applies  the  shifts  from  a  refinement  cycle  to 
the  input  parameter  file  (unit  10),  with  damping  factors, 
if  so  specified,  and  writes  the  new  parameter  file  onto 
unit  11.  In  addition,  it  calculates  and  prints  the  average 
position  shift  (in  Angstroms)  before  and  after  damping, 
keeps  track  of  the  largest  atomic  shift  and  prints  it  out, 
and  currently  flags  all  atoms  which  move  more  than  0.25  A 

with  quadruple  astericks.  Also,  if  the  scale  factor  changes 
by  more  than  20%  of  itself,  the  thermal  factors  are  not  refined 

at  all  that  cycle.  (There  is  a  strong  correlation  between  scale 
factor  and  the  overall  magnitude  of  the  thermal  factors,  which 
can  lead  to  lengthy  oscillations  in  the  refinement.  Damping  of 
the  scale  and  thermal  factor  shifts  is  probably  always  advisable 
we  usually  use  a  damping  factor  of  0.5  ^»ven  near  convergence). 

8.  RZSLSQ  is  the  main  program,  the  central  driving  program. 

It  doesn't  really  do  much  except  sequentially  call  each  of  the 
major  subroutines,  and  print  out  a  summary  of  times  once  each 
cycle.  It  might  be  useful  for  getting  an  idea  of  the  overall 
program  flow. 

Note:  Its  last  call  is  to  RECYL  -  this  is  an  alternate 

entry  point  of  INIT. 

9.  DI AG  -  Estimates  the  standard  deviations  for  the  parameters 
being  varied.  The  S.D.'s  are  based  only  on  the  diagonal  matrix 
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elements  (the  contributions  from  the  restraints  are  subtracted 
out)  . 

10.  DISRF  calculates  the  contirbution  of  the  distance  restraints 
to  the  residual  and  to  the  normal  matrix  coefficients.  Prints 
out  a  summary  of  actual  and  target  distances,  with  large 
discrepancies  flagged.  This  should  be  checked  carefully  in  the 
beginning  to  guard  against  inadvertent  connectivity  errors.  In 
later  cycles,  watch  the  actual  SIGD  versus  the  specified  SIGD  for 
the  various  classes.  Ideally,  they  should  roughly  match.  If 
actual  sigmas  are  several  times  larger  than  specified  sigmas,  the 
x-ray  data  may  be  trying  to  tell  you  your  target  distances  are 
not  so  good,  and,  in  fact,  the  targets  may  be  keeping  you  from  a 
better  R-factor.  REDIST  is  an  entry  into  DISRF  which  is  used  to 
skip  unnecessary  calculations  when  only  an  R-factor  is  wanted. 

11.  ELAPSD(N)  -  A  call  to  ELAPSD  will  zero  all  eight  time 
counters,  TIME(l)  to  TIME(8),  if  N=0.  (It  should  be  initialized 
in  this  way,  and  is,  in  RESLSO.)  Later,  when  ELAPSD  is  called, 
the  elapsed  time  since  the  last  previous  call  is  added  to  TIME (N) 
This  is  accomplished  on  our  computer  by  a  special  system  function 
S$SET(SCCB),  which  returns  "elapsed  central  processor  time  in 
seconds",  each  time  it's  called. 

12.  FLAG  (I,N,A,B,C)  -  A,  ideal,  is  compared  with  3,  real,  and 
divided  by  C,  the  est.  stand,  dev.  of  diff.  N  returns  a  blank,  a 
single  asterisk,  or  a  triple  asterisk  depending  on  discrepancy 
magnitude . 
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13.  FREEZE  -  If  the  HRIDE  option  has  been  chosen,  then  FREEZE 
will  set  the  matrix  and  vector  elements  for  hydrogen  atom  thermal 
factors  to  zero.  This  essentially  removes  the  hydrogen  thermal 
parameters  from  the  refinement. 

14.  HRIDE  -  This  subroutine  is  only  used  if  IHRID=1  which 
indicates  that  hydrogen  thermal  parameters  are  not  to  be  refined. 
HRIDE  searches  the  atoms  list  for  hydrogens  and  then  identifies 
the  atoms  to  which  they  are  bonded.  It  then  sets  the  orientation 
of  the  thermal  ellipsoids  for  the  hydrogens  to  be  the  same  as  those 
of  the  heavier  atoms  to  which  they  bonded  and  it  sets  the  ellipsoid 
size  to  be  equal  to  or  greater  than  (by  percentage  BRAT)  those  of 
the  heavier  atoms. 

15.  INIT  -  We  have  tried  to  localize  all  reading  of  card  input 
in  INIT  (File  26  is  read  in  DISRF) .  Also,  most  initializing 

of  arrays  and  tables  is  done  here.  Entry  RECYCL  (near  the  end) 
rewinds  and  re-initializes  units  and  arrays  in  preparation  for 
further  refinement.  ' A  recap  of  all  the  refinement  control  para¬ 
meters  is  printed  out.  The  comment  cards  embedded  in  this 
subroutine  explain  many  of  the  possible  refinement  options  available 
to  the  user. 

15  &  17  -  MATMLL  and  MATMULL  -  These  routines  do  a  simple 
multiplication  of  a  matrix,  A,  times  vector,  P,  to  yield  a 
vector,  F.  However,  it  is  not  so  simple  because  A  is  not  stored 


as  a  2-dimensional  array  (because  it  is  sparse,  etc...).  To 
do  the  multiplication,  a  "bookkeeping  scheme"  must  be  used  to 
determine  which  elements  of  A  multiply  the  components  of  P. 

MATMUL  is  used  for  the  position  matrix,  MATMLL  for  the  thermal 
factor  matrix.  One  would  suffice,  with  some  re-coding. 

18.  PLANE  uses  the  method  of  Schomaker,  et  al.  (Acta  Cryst. 
(1959)  12,  600)  to  determine  the  best  least-squares  plane  through 
a  set  of  points. 

19.  PLNRF  finds  the  deviations  of  atoms,  in  groups  specified 
as  planar,  from  a  best  least-squares  plane.  It  calculates  and 
prints  out  the  individual  deviations,  and  their  contributions  to 
the  residual,  and  increments  the  normal  matrix  coefficients  and 
vector  elements. 

20.  RECYCL  (see  INIT) 

21.  REDIST  (see  DISRF) 

22.  RFAC  -  A  "summary"  routine.  Subroutine  SF  forms  12  sums, 
which  are  (misleadingly)  labelled  SIG(l)  to  SIG(12) .  By  com¬ 
bining  these  12  sums  in  various  ways  to  form  numerators  and 
denominators,  an  incredible  array  of  R  factors  are  formed,  here 
in  RFAC,  and  printed  out.  Fairly  straight  forward. 

23.  RSCAL  (see  SCALE) 

24.  SCALE  (see  also  CG30LV)  -  SCALE  is  used  to  perform  a  trans¬ 
formation  to  the  normal  equations: 


where  S  is  a  diagonal  scaling  matrix,  and  its  components  are 
adjusted  so  as  to  equalize  the  diagonal  components  of 
A'  =  SAS 

The  conjugate  gradient  procedure  is  used  to  solve  the  new 

set: 

A 1 x '  =  r  where  x 

and  ?  =  Sr 

Finally,  the  desired  shift,  x,  is  retrieved  by  calling 
entry  RSCAL  and  calculating  the  product: 

X  =  Sx' 

By  equalizing  the  diagonal  components,  SCALE  transforms  to 
a  set  of  normal  equations  in  which  the  derivative  of  the  parameters 
employed  in  the  Taylor  expansion  are  re-weighted  so  that  they 
have  equal  average  values.  Thus  equal  components  in  the  S~^  vector 
have  equal  structural  significance.  This  aids  in  the  convergence 
of  the  conjugate  gradient  solvent  to  the  normal  equation,  since 
initial  iterations  determine  the  largest  components  of  S  most 
accurately. 

25.  SCALEB  -  Same  operations  as  SCALE  but  for  thermal  shifts. 

26.  SF  -  Subroutine  SF  reads  the  reflection  data  from  unit  IREFR. 
It  checks  to  see  whether  a  reflection  is  to  be  included  in  the 
refinement.  If  it  is,  CALC  is  called  to  compute  FCALC.  Weights 
are  applied  here  and  the  various  agreement  sums  are  incremented. 

SF  calculates  phases  and  writes  the  difference  map  (or  structure 
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factor  listing)  file  as  desired.  SF  also  does  the  derivative 
calculations  and  increments  the  matrix  and  vector  elements  in  their 
respective  arrays.  (This  part  of  SF  is  skipped  when  only  an  R- 
factor  is  desired.) 

27.  SYMMIN  -  The  main  function  of  this  subroutine  is  the  calcula¬ 
tion  of  some  Cartesian-based  symmetry  operation  rotation  matrices, 
labelled  ARA  (3, 3, NS) ,  which  later  are  used  in  CALC  to  determine 
the  ellipsoid  orientation  matrices  for  symmetry  mates  of  the  basic 
atom  deck. 

23.  TTST  -  Checks  to  see  if  any  thermal  factors  have  gone  negative 
and,  if  so,  flags  them  but  does  not  stop  the  refinement. 


RESLSQ  INPUT  FILES 


FILE  05 
FILE  10 


■C  i. 
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CONTROL  CARDS  (described  elsewhere) 

ATOMIC  COORDINATES  (read  in  SUBROUTINE  INIT) 

1. )  scale  factor  card  FORMAT  (F10.5) 

2. )  coordinate  cards  -  2  per  atom 

Card  1:  FORMAT  ( A8 , IX, 5F9 . 6 ) 


ATOM 

FMULT 

OC 

XYZ  (1) 
XYZ (2) 
XYZ  (3) 


atom  name 

scattering  factor  reference 

occupancy 

x 

Y 

2 


Card  2:  FORMAT  (6F9.6) 


BETA ( 1) 
BETA (2) 
BETA ( 3 ) 
BETA (4) 
BETA (5) 
BETA (6) 


thermal  parameters 


REFLECTION  AND  SCATTERING  FACTOR  DATA 
(Read  in  SUBROUTINE  SF) 


BINARY  FILE  containing 


IH1  -  h 
IH2  -  k 
IH3  -  I 

YO  -  F  OBSERVED 
SIGG  -  WEIGHT 
STHOL  -  SIN  0/X 
FII (NNG) -  scattering  factor 

(one  for  each  type  of  atom  in  the  molecule) 


NOTE : 


Each  reflection  carries  its  own  scattering  factor  informa¬ 
tion.  This  helps  speed  up  structure  factor  calculations. 
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FILS  26  DISTANCE  RESTRAINT  PARAMETERS  (read  in  SUBROUTINE 

DISRF)  3INARY  FILE 

NOR  -  ATOM  NUMBER  OF  ORIGIN  ATOM 
NT  -  ATOM  NUMBER  OF  TARGET  ATOM 
DIS  -  IDEAL  DISTANCE  BETWEEN  NOR  AND  NT 
NWTD  -  BOND  TYPE 

1  =  30NDED  DISTANCE  (NON-H) 

2  =  NON-BONDED  DISTANCE  (NON-H) 

3  =  BONDED  DISTANCE  (INVOLVING  H) 

4  =  NON-BONDED  DISTANCE  (INVOLVING  H) 

5  =  INTERMOLECULAR  DISTANCE 


XSYMMl 

XSYMM2  LATTICE  TRANSLATIONS  FOR  TARGET  ATOM 
XSYMM3 

LSYMl  SYMMETRY  OPERATOR  NUMBER  FOR  TARGET  ATOM 
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RESLSQ  -  OUTPUT  FILES 


i 


FILE  06  PRINT  OUT 

FILE  07  (Optional;  written  by  BCNVRT) 

Output  file  of  standard  deviations  for  coordinates 

1  card  per  atom;  FORMAT  (A8,  19X,  3F9.6) 

atom  name 

SDX 

SDY 

SDZ 

FILE  11  (written  by  COMBB) 

Output  file  for  shifted  coordinates.  Scale  factor 
is  written  on  the  first  card,  then  two  cards  per 
atom  -  same  format  and  variables  as  INPUT  FILE  10 

FILE  30  (written  by  SF)  FORMAT  (1X,3I4  ,  5F11 . 4 , 11X, Fll . 4 , 

7X, F10 .5) 

Structure  factor  file 

IH  ( 1)  h 

IH  (  2 )  k 

IH ( 3)  1 

YO  FOBS 

YC  FCALC 

APHAS  A 

BPHAS  B 

DY  FOBS -FCALC 

SIGYOA  WEIGHTED  (FOBS-FCALC) 

XF  PHASE  (IN  RADIANS) 
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The  ASC  computer  and  possible  machine-specific  code 


I.  Functions  used 


I  =  ..’D  (L,M) 


£(I  =  L.  AND  . M) 
I  =  MINO  (L,M) 


AMIN1  (A, B) 


Performs  a  "Boolean  and"  operation  on  L  and 


M. 

for  VAX  version  of  the  program 


•3 


Finds  lowest  of  integer  arguments  (L,M)  and 
places  it  in  I. 

Finds  lowest  of  floating  point  arguments. 


The  ASC  has  a  32  bit  word.  I  don't  think  this  is  highly 
pertinent  to  the  operation  of  the  code  in  RESLSQ,  except  that 
roundoff  error  is  appreciable,  and  one  shouldn't  expect  to 
reproduce  test  outputs  exactly. 

The  ASC  has  a  vectorizing  capability  -  a  special  compiler 
and  hardware  for  processing  calculation  on  arrays  rapidly.  This 
has  affected  (sometimes  more  than  necessary)  the  way  we  set  up 
our  calculations  in  FORTRAN,  Simple  DO-loops  without  IF-checks 
will  generally  vectorize  (and  if  the  loop  is  fairly  long,  this 
can  sometimes  increase  speed  by  a  factor  of  2C  or  30!).  Thus, 
the  coding  may  sometimes  seem  curiously  strung  out  into  many 
successive  DO-loops,  when  it  might  be  more  compact.  This 
shouldn't  hamper  its  operation,  and  probably  not  its  speed,  on 
other  compilers  or  computers. 

ASC  memory  is  large  and  cheap  at  present  -  sc  the  program 
is  lavishly  dimensioned.  3esides  lowering  dimensions,  it  might 
be  possible  to  do  away  with  small  arrays  through  re-coding  (See 
SCALE) .  Also,  more  overlaying  of  space  with  EQUIVALENCE  statements 
may  be  possible. 
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Implied  DO-loops  -  In  ASC  FORTRAN,  one  can,  within 
certain  limits,  use  array  names  without  explicit  subscripts 
in  arithmetic  statements  and  the  statement  will  be  interpreted 
as  applying  to  the  whole  array. 

Example 

DIMENSION  A(900) ,  B(900) 

B  =  0.0 
A  =  B  +  2.0 

The  first  statement  is  equivalent  to 
DO  10  J  =  1,900 
10  B ( J)  =  0.0 

and  the  second 
DO  20  J  =  1,900 
20  A ( J)  =  B ( J)  +  2.0 

We  have  tried  to  avoid  using  special  FORTRAN;  however,  the 
above  feature  is  so  convenient  that  it  might  possibly  occur  some 
where.  If  so,  the  effects  may  be  treacherous.  We  would  like  to 
know  if  any  instances  are  discovered. 


ARRAYS  IN  SUBROUTINE  DIMENSIONS 


BRF : 


BCNVRT : 


New  Locations 


Til (3,6) 

18 

AOR (3,6) 

18 

AP (3,6) 

18 

DP (3,3) 

9 

DIR  ( 8 ) 

8 

NSUM ( 8 ) 

8 

RMSD (8) 

8 

Bll (500) 

EQUIVALENCE 

AM  (1) 

B22 ( 500 ) 

EQUIVALENCE 

AM  (501) 

B33  (500) 

EQUIVALENCE 

AM (1001) 

B12 (500) 

EQUIVALENCE 

AM (1501) 

B13  (500) 

EQUIVALENCE 

AM  (2001) 

B23 (500) 

EQUIVALENCE 

AM (2501) 

BUS D  (500) 

EQUIVALENCE 

AM (3001) 

B22SD  (500) 

EQUIVALENCE 

AM  (3501) 

B33SD  (500) 

EQUIVALENCE 

AM (4001) 

B12SD  (500) 

EQUIVALENCE 

•AM  (4501) 

B13SD  (300) 

EQUIVALENCE 

AM (5001) 

B23SD  (500) 

EQUIVALENCE 

AM (5501) 

IDX  (3,500) 

EQUIVALENCE 

AM  (6001) 

BEQSD  (500) 

EQUIVALENCE 

AM (7501) 
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New  Locations 


CALC : 

TT (3,3,500) 

9 

X  NA 

ARG ( 500 ) 

NA 

DFDS  (18) 

SC 

COS1  (500) 

NA 

SIN1  (500) 

NA 

FI  (500) 

NA 

INDX  ( 500 ) 

NA 

DA (3,500) 

3 

x  NA 

TFl  (500) 

NA 

TFX1  (500) 

NA 

HH  (30) 

3 

19xNA  +  SC  +  3 

CART : 

XX (3,500) 

3 

x  NA 

X ( 3 , 500 ) 

3 

x  NA 

CART ( XX , X , NTC ) 

CGSOLV : 

P (1518) 

3 

x  NA+SC 

EQUIV (DADX (1,1) ) 

F (1513) 

3 

X  NA+SC 

EQUI7 ( DBDX ( 1 , ) ) ) 

CGS0L1 : 

P( 3001) 

6 

x  NA+SC 

EQUIV (DADX  (1, 1) ) 

F( 3001 ) 

6 

x  NA+SC 

EQUIV (DBDX (1, 1) ) 

COMBB : 

XYZN  (3, 5  00) 

x  NA 

EQUIV (DADX (1,1) ) 

BETAN (3,500) 

3 

x  NA 

EQUIV (DAD3 (1,1)) 

1FLAG ( 500 ) 

NA 

EQUIV (DBDX (1,1)  ) 

BEQ ( 500 ) 

NA 

EQUIV (DBDB  (1,1)) 

PLNRF : 


3  x  NA  +  SC 


6  +  4  x  NAP 


9 

+  SC 

+  SC  SCALE (A ,  FK) 

+  SC 


BK( 1518) 
DEL (20) 
VM(B  ) 

X( 3 . 20 ) 
DIR( 3 ) 

RFAC:  RCH ( 4 ) 

SUM( 5 ) 

SCALE:  FK( 1518 ) 

X( 1518) 
SCL( 1518 ) 
*A( 40000 ) 
*S  ( 40000 ) 

SCALE  B:  FK(3001) 

X( 3001 ) 
SCL( 3001 ) 
*A( 40000) 

*S ( 40000 ) 


NAP  PLNRF (BK) 

3 

3  x  NAP 

3 

4 

5 

3  x  NA 
3  x  NA 
3  x  NA 
40000 
40000 

6  x  NA  +  SC 
6  x  NA  +  SC 
6  x  NA  +  SC 


40000  +  6  x  NA  +  6  X 


40000  +  12  x  NA  +  12  X  SC 


NOTE:  THE  S  ARRAY  HERE  IS  USED  SOLELY  TO  ENABLE  VECTORIZING 

CERTAIN  CALCULATIONS  IN  SCALE  A  AND  B.  IT  COULD  BE 
RE-CODED  TO  ELIMINATE  THE  40000  LOCATIONS  USED  UP  BY  S. 
IN  THE  VAX  VERSION  OF  THE  PROGRAM  THIS  ARRAY  HAS  BEEN 
ELIMINATED. 


SC 


*Not  big  enough  for  500  atoms 


•'  -*■  "  -*'*''*  -*  ^ 

11 

SF: 

INDX( 2500 ) 

NDIS 

DVO( 2500) 

NDIS 

"  V 

DVZ( 2500 ) 

NDIS 

I 

IHKL( 4 ) 

4 

NSG( 10) 

10 

14  +  3  x  NDIS 

CONLSG : 

NO  DIMENSION 

STATEMENTS 

D I  AG : 

NO  DIMENSION 

STATEMENTS 

DISRF : 

DSUM( 8 ) 

8 

NSUM( 8 ) 

8 

RMSD ( 8  ) 

8 

a 

DIR(6) 

6 

30 

ELAPSD : 

SCCB( 3 ) 

s 

i  , 

3 

FLAG : 

NO  DIMENSION 

STATEMENTS 

INIT: 

SPACGRP( 5 ) 

5 

•# 

* 

MTAMLL : 

P( 1518 ) 

3  x  NA  + 

SC 

MATMLL ( P , F ) 

5 

1  • 

F ( 1518) 

3  x  NA  + 

SC 

- 

MATMUL : 

P( 1518) 

3  x  NA  + 

sc 

o 

F  ( 1 5 1 8 ) 

3  x  NA  + 

sc 

MATMUL ( P , F ) 

58 

— 

0 

# 

PLANE : 


9 


G(3,3) 

X( 3 , 20)  60 

XS( 3  )  3 

XXS (3,3)  9 

ADJ (3,3)  9 

B( 3 , 3  )  9 

A( 3 , 3 )  9 

VM ( 3 )  3 

VMI ( 3 )  3 

BV ( 3 )  3 

SYMMIN :  TEMP (3, 3)  9 

AS(3 , 3)  9 

SI ( 3 , 3 )  9 


UNIT:  XXX( 3 , 500)  3  x  NA 

DIST( 500 )  NA 

TOTAL  LOCATIONS  CURRENTLY  USED  BY  DIMENSIONS  STATEMENTS 


41  x 

NA  = 

20500 

19  x 

SC  = 

95 

4  x 

NAP  = 

120 

3  x 

NDIS  = 

7500 

MI  SC 

80309 

108524 
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